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Abstract 

This paper investigates the impact of parameter uncertainty on capital estimate 
in the well-known extended Loss Given Default (LGD) model with systematic 
dependence between default and recovery. We demonstrate how the uncertainty 
can be quantified using the full posterior distribution of model parameters ob- 
tained from Bayesian inference via Markov chain Monte Carlo (MCMC). Results 
show that the parameter uncertainty and its impact on capital can be very sig- 
nificant. We have also quantified the effect of diversification for a finite number 
of borrowers in comparison with the infinitely granular portfolio. 
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1 Introduction 



Default and recovery rates are key components of Loss Given Default (LGD) models 
used in practice for calculation of economical capital (EC) against credit risk in some 
major banks. The classic LGD model implicitly assumes that the default rates and 
recovery rates are independent (Bluhm et al 2002). Frye (2000), Pykhtin (2003) and 
Dullmann and Trapp (2004) extended the classic model to include systematic risk in 
recovery rates, incorporating a non-zero correlation between default rates and recovery 
rates driven by the systematic factor. They considered three extensions to account 
for the systematic risk in recovery rates under three different assumptions for the 
distribution of recovery rates: Frye (2000) - a normal distribution; Pykhtin (2003) 
- a log-normal distribution; Dullmann and Trapp (2004) and Schonbucher (2001) - 
a logit-normal distribution for the the recovery rate. The extended models are still 
parsimonious, yet they represent an important enhancement of credit risk models used 
in earlier practice, for example, CreditMetrics (Gupton et al 1997) and CreditRisk+ 
(Credit Suisse Financial Products 1997) that do not account for systematic risk factor 
driving both default and recovery rates. 

Dullmann and Trapp (2004) summarized the empirical literature on systematic risk 
in recovery rates, and found there was a broad agreement that default rates and business 
cycle are correlated. They provided closed-form solutions for the maximum likelihood 
estimators (MLEs) of parameters for the probability of default and recovery rate distri- 
butions. They also estimated the correlations of default rates and recovery rates with 
the systematic risk factor, and demonstrated the consequences for EC if the classic 
one-factor model is applied instead of their extended models. In particular, it was 
observed that EC is significantly higher in the extended models due to dependence of 
recoveries on the systematic risk factor; and EC estimates are close when comparing 
results under three different assumptions for the recovery rate distribution mentioned 
in the previous paragraph. 

Typically, the available data consist of default and recovery rates for the time period 
concerned (e.g. one year) and are of limited size, covering a couple of decades at most. 
For example, in the study of Dullmann and Trapp (2004), the default and recovery data 
have eighteen points covering an eighteen-year period. Inevitably the limited data 
size could pose significant instability and uncertainty in the LGD model parameter 
estimations. None of the various studies, including the extension work of Frye (2000), 
Pykhtin (2003) and Dullmann and Trapp (2004) specifically addressed the quantitative 
impact of parameter uncertainty. 

Increasingly, quantification of parameter uncertainty and its impact on capital has 
become a key component of financial risk modeling and management. Recent examples 
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of addressing parameter uncertainty in operational risk and insurance capital modeling 
include Luo et al (2007) and Peters et al (2009a). This paper studies parameter un- 
certainty and its impact on EC estimate in the extended LGD model, i.e. probability 
of default model coupled with a dependent recovery rate model. We demonstrate how 
the model with latent systematic factor can be estimated using Bayesian approach and 
Markov chain Monte Carlo (MCMC) method. 

The organization of this paper is as follows. Section 2 first describes the credit 
risk model setup, particularly the extended default and recovery models considered by 
Frye (2000), Pykhtin (2003) and Diillmann and Trapp (2004). This is followed by a 
discussion on various EC estimates and the corresponding algorithms, either for an 
arbitrary portfolio and for the limiting case of a infinitely granular portfolio. Section 3 
presents the probability and likelihood functions for the model. This includes the full 
joint likelihood for default and recovery; as well as the approximate likelihood method 
used in Diillmann and Trapp (2004). Section 4 describes the Bayesian inference formu- 
lation and the MCMC simulation algorithm for the posterior distribution of the LGD 
model parameters. Section 5 presents MCMC results in comparison with the MLEs, 
demonstrating the impact of parameter uncertainty on economic capital estimation. 
The effect of diversification is also shown in Section 5. Concluding remarks are given 
in Section 6. 



The classic one-factor model assumes a homogenous loan portfolio where the distribu- 
tion of its loss vector that collects losses of individual loans is exchangeable, or invari- 
ant under permutations of its components. Following Frye (2000), Pykhtin (2003) and 
Diillmann and Trapp (2004), the key characteristics of the classic one-factor model are 
summarized as follows. 

Consider a portfolio of J borrowers (firms) over a chosen time horizon, where the jth 
borrower has one loan with principal amount A,-. The loss rate (loss amount relative 
to the loan amount) of the portfolio due to defaults is 



2 LGD Model 




(i) 



where we have the following definitions: 



• Wj - the weight of loan j in the portfolio, Wj = 




• Lj - the loss rate of loan j due to default; 
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• Rj - the recovery rate of loan j after default; 

• Ij - indicator variable associated with the default of firm j; Ij = 1 if firm j 
defaults, otherwise Ij = 0. 

Formally, Rj should have [0, 1] support. However, some distributions proposed for 
Rj have support beyond [0,1], e.g. the normal distribution. Indeed, in practice, the 
recovery may appear beyond [0, 1] interval. Following Diillmann and Trapp (2004), in 
this study we do not explicitly impose the restriction < Rj < 1. In fact, results in 
Diillmann and Trapp (2004) show that the unbounded normal distribution for recovery 
rate gives a capital estimate very close to that given by the properly bounded logit- 
normal distribution (the relative difference is less than 1%). 

Remark: The above notation is for a given time period. Later, starting from Section 
3, we consider the model over a number of time periods t — 1, 2, . . . , T, T + 1 that will 
add index t to all random variables. Here, T + 1 refers to the next year. It is assumed 
that all random variables involved in the model are independent between different time 
periods. However, the model can be easily extended to have explicit time dependence. 

2.1 Modeling Default 

Denote the probability of default for firm j by p, i.e. Pr[ij = 1] = p. Let Cj be an 
underlying latent random variable such that firm j defaults if Cj < where 
$(•) is the standard normal distribution and < I )_1 (-) is its inverse. That is, Ij = 1 
if Cj < and Ij = otherwise. Cj describes the overall financial condition 

(financial well-being) of firm j over a time horizon. The value Cj for each firm depends 
on a systematic risk factor, denoted by X, and a firm specific (idiosyncratic) risk factor 
Zfas 

Cj = ^-pX+ y/T=jZf, (2) 

where 2Tf , . . . , Zj are pairwise uncorrected. Also X and Z^f are assumed to be inde- 
pendent and from the standard normal distribution. 

Conditional on X, the financial conditions of any two firms are independent. The 
parameter p quantifies the extent of exposure of a firm's asset value to the fluctuations 
in the business cycle. Unconditionally it measures the correlation between financial 
conditions of two firms. The value of p, < p < 1, is assumed to be constant for all 
firms at all time periods. 
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2.2 Recovery following default 

The extended models considered by Frye (2000), Pykhtin (2003) and Diillmann and 
Trapp (2004) account for systematic risk in recovery rates under three different as- 
sumptions for the distribution of recovery rates. Define 



Vj = fJL + CTy/ujX + oV 1 - tuZj, [0,1], (3) 

where X and Zj are assumed to be independent from the standard normal distribution, 
and the parameter u is restricted to the interval [0, 1]. Also, Zj and Zj 7 are assumed 
independent too. Note, the one-factor model in Pykhtin (2003) allows for correlation 
between Zj and Z^ . The three models for the recovery rate are then defined through 
Vj as follows. 

• The first extended model, as initially suggested by Frye (2000), assumes a normal 
distribution for the recovery rates, i.e. the recovery rate Rj of obligor j is given 
by 

Rj = Vj. (4) 

An advantage of the above model is that the parameters fi and u directly represent 
the mean recovery and recovery correlation. 

• The second extension, initially proposed by Schonbucher (2001), assumes that 
the recovery rate Rj follows a logit-normal distribution, i.e. 

Rj = eMv j, . (5) 

The above model satisfies the restriction < Rj < 1. 

• The third model, following Pykhtin (2003), has a log- normal distribution for the 

recovery rate 

Rj = exp(Vj). (6) 

One of the the findings by Diillmann and Trapp (2004) is that the EC estimate is 
insensitive to the model choice for the above recovery models - only about 2% difference 
was observed among the EC values estimated by these three models. In addition, they 
carried out Shapiro- Wilk test and Jarque-Bera test for normality, and found that the 
normal distribution assumption for the recovery rate is favored by the p-values over 
the other two models. Thus in the present study, we will concentrate on the first 
recovery model given by (jl]), i.e. we assume a normal distribution for the recover rate. 
Another reason for our choice of model (HI) is because we do not have the data for 



5 



individual recoveries but only the average recovery rates; and we can use the fact that 
the distribution of the average of normally distributed independent random variables 
is still normal. 

2.3 Economic capital 

The EC for credit risk can be estimated using a high (e.g. 0.999) quantile for the 
distribution of loss L defined in (1). Specifically, the quantile Q q is defined as 

Q g (0) =Q q = mi{z : Pr[L > z] < 1 - q} = mi{z : F L {z) > q}. (7) 

where q is a quantile level (e.g. 0.999) and Fl(z) is distribution function for the random 
loss L; the corresponding density of L is denoted as The dependence of Fl(z) 

on the model parameters 6 = (p, p, fi, a, u>) is implicit in (J7J). 

There are different ways of estimating this high quantile, some are based on point 
estimates (e.g. MLEs) of parameters and others account for parameter uncertainty. We 
are interested in comparing these different estimates of capital charge and quantifying 
the impact of different assumptions, particularly the impact of parameter uncertainty. 

2.3.1 Point estimates of capital 

For a given model with parameters 6, the EC defined as a quantile (jZJ) is a function of 
6. Typically, given observations, the MLEs 6 are used as the "best fit" point estimators 
for 0. Then, the loss density for the next time period is estimated as f L (z\6) and its 
quantile, Q q (6) } is used for the economic capital calculation. In general, the distribution 
of L is not tractable in closed form for an arbitrary portfolio. In this case Monte Carlo 
simulation can be used for the capital estimation. Briefly, the procedure for simulating 
L in (1) with given parameters can be described as follows. 

Algorithm 1 (Quantile given parameters) 

1. Draw a single independent sample from the standard normal for the systematic 
factor X. 

2. For each borrower j = 1, . . . , J: draw an independent sample from the standard 
normal distribution for the idiosyncratic default risk factor Zj, compute Cj as 
in (2); and let Ij = 1 if Cj < $ _1 (p) and Ij = otherwise. 

3. Draw an independent sample from the standard normal distribution for the id- 
iosyncratic recovery factor Zj and compute Rj — o^fojX + ay/1 — uZj. 
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4. Compute loss L for the entire portfolio using (1). This is a sample from the loss 
distribution. 

5. Repeat steps 1-4 to obtain iV samples of L with iV sufficiently large for high 
quantile calculations. 

6. Estimate Q q {6) using obtained samples of L in the standard way. 

In practice, the parameters are unknown and it is important to account for this 
uncertainty when the capital charge is estimated, especially for small datasets. A 
standard frequentist approach to estimate this uncertainty is based on limiting results 
of normally distributed MLEs for large datasets. Then information matrix is used to 
estimate the covariances between MLEs. In this paper we take Bayesian approach, 
because dataset is small and the distribution of parameter uncertainty is very different 
form normal. Estimation of the quantile accounting for parameter uncertainty under 
a Bayesian framework will be discussed in Section 4. 

2.3.2 Economic capital under the limiting condition 

In the case of a diversified portfolio with a large number of borrowers, the idiosyncratic 
risk can be eliminated and the loss depends on X only. Gordy (2002) has shown that 
the distribution of portfolio loss L has a limiting form as J — > oo, provided that each 
Wj goes to zero faster than 1 / \f~J . The limiting loss rate L°° is given by the expected 
loss rate conditional on the systematic factor X 



j j 

L°° = L°°(X) = E[L\X] = ^wjEfalX] = ^V^max(l-i^0)|X], (8) 

3=1 3=1 

i.e. the limiting loss L°° is just a function of X and the distribution of L°° is fully 
implied by the distribution of X. 

Conditional on X the default indicator variable Ij and the recovery rate Rj are 
independent, since in (2) and Zj in (3) are independent. Thus limiting loss flH]) for 
J — > oo becomes 

j j 
L°° = ^w j E[I j \X}E{m^(l - Rj,0)\X] = J2 w j A j( x ) S j( x )> ( 9 ) 

3=1 3=1 

where Aj(X) = E[Ij\X] is the conditional probability of default of firm j and Sj(X) = 
£ , [max(l — Rj,0)\X] is the conditional expected value of loss rate, both are functions 
of X. 
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Bank loans are subject to borrower-specific risk and systematic risk. The former 
can be controlled or even neutralized by diversification. For a homogeneous portfolio, 
probability of default and recovery rates (or loss given default) are not firm specific, 
i.e. Aj(X) = A(X) and Sj(X) = S(X) for all j, and © simplifies to 

j 

L°°(X) = J2 w j A (X)S(X) = A(X)S(X). (10) 

i=i 

That is, the limiting loss rate of the diversified homogenous portfolio is a function of 
X only. As in the model underlying the internal ratings-based risk weights of Basel II, 
EC is determined with the assumption that the bank loan portfolio is fully diversified 
and the EC is only held for systematic credit risk. It can be readily shown that the 
quantile of L°°(X) at level a can be calculated as Q™ = L°°(X = $ _1 (1 — a)) because 
L°°(X) is a monotonic decreasing function of random variable X and X is from the 
standard normal distribution. As in Dullmann and Trapp (2004), we define the EC of 
the diversified portfolio loss distribution L°°(X): 



EC™ = Q °° 999 = ^(^(O.OOl)) 

= A($ _1 (0.001)) x ^((^(O.OOl) = PD x LGD, (11) 

where PD = A($- 1 (0.001)) and LGD = S($- x (0.001). Given X, the conditional 
probability of default, A(X), can be found from (2) as 



A(X) = $ 



:i2i 



The expected conditional loss rate for the normally distributed recovery rate model 
(4) is calculated in Dullmann and Trapp (2004) as 

S(X) = £[max(l -Rj,0)\X] « E[(l - Rj)\X] = 1 - fi - a^/uX. (13) 

The approximation E7[max(l — Rj, 0)\X] ~ E[(l — Rj)\X] assumes that Rj exceeding 1 
does not have much impact on the results. At least for the set of real data used in the 
present study, this approximation makes no material difference. It is not difficult to 
avoid this approximation either here or in any of the Monte Carlo simulations. In order 
to compare with results of Dullmann and Trapp (2004) under the same assumptions, 
we use formula ( ITBl . As noted before, in this study the conditional probability of 
default A(<3> -1 (0.001)) is denoted as PD, and the corresponding conditional loss rate 
^(^(O.OOl)) is denoted as LGD. 
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Under the framework outlined above, the EC°° is a function of the five parameters 
= (p, p, //, a, co), with X = $ _1 (0.001) ~ —3.09. Obviously an uncertainty in any of 
the the parameter estimates will cause an uncertainty in the EC estimate. This will 
be discussed in Section 4.1. 

Dullmann and Trapp (2004) calculate the economical capital using ( II ip as EC°° = 
Qo D Qgg(0 MLE ), where Q MLE is the point estimate of parameters obtained via approxi- 
mate maximum likelihood method described in Section 13.21 

3 Likelihood 

3.1 Exact likelihood function 

Consider time periods t = 1,2, ... ,T, T+l with T+l corresponding to the next (future) 
time period, where the following data of default and recovery for a loan portfolio of J t 
firms are observed: 

• D t - the number of defaults in time period t, with d t denoting the actual realiza- 
tion observed; 

• \l/t - the default rate in time period t, = D t / J t , with ip t denoting the actual 
realization observed; 

• Rt - the average recovery rate in time period t, with f t denoting the actual 
realization observed. 

Denoting the individual recovery rates for the d t defaulted firms as R\{t), . . . , R<i t {t), 
the observed average recovery rate is R t = J2f=i Rj{t)/dt- The joint density of {D t , R t ) 
can be computed as 

f(d t ,r t ) = J f(d t \x t )f(f t \d t ,x t )f N (x t )dx t . (14) 

Given X t = Xt, all the firms in a homogenous loan portfolio have the same conditional 
default probability Pr[Ij(t) = l\X t = x t ] = A(x t ) evaluated in ( ITS]) . Since D t = 
Y2jLi Ij(t)> the distribution of D t is a binomial, that is 

f(d t \x t ) = Pr[A = dt\X t = x t ] = Q) (A(x t )) dt (1 - A(x t )) Jt - dt , (15) 

which can be well approximated by the normal distribution N(/j, t , of) for large J t , 
where \i t = JtA(x t ) and of = J t A(x t )(l — A(x t )). That is 
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Conditional on X t — x t , Rj(t) (j = 1, . . . , d t ) are independent from normal distribution 
N(p r , of) with \x r = fi+a^/uxt and ay = cry/1 — uj. Thus the average R t is from normal 
distribution N(fi R , cr R ) with fi R = \i r and a\ = cr^/d t 

f(f t \D t = d u X t = x t ) = -j^exp f- (rt ~f H . (17) 

Substituting f|T6|) and ( 1T7|) into f lT4|) . the joint density f(d t ,r t ) can be calculated nu- 
merically. Define random vectors D = (Di, . . . , D?) and R = (Ri, . . . , -Rt)- The joint 
likelihood function for D and R is then 

T 

^ D ^0) = Hf(d t ,r t ). (18) 
t=i 

The likelihood function fl 18p involves numerical integration in (114p . which integrates 
out the latent variables X±, . . . ,X? and thus reducing the dimension (or number of 
unknown parameters) of the problem significantly (from 23 to 5 for the dataset con- 
sidered in this study). However, our numerical experiments show that, although not 
impossible, it is difficult in practice to accurately compute integrations in ( 1T8|) . One of 
the difficulties is the frequent occurrence of numerical under-flow in the evaluation of 
the integrand, even in double precision using Gauss-Hermite quadrature. 

A more straightforward and problem-free approach (avoiding integration) is to treat 
the latent variable X the same as other parameters, and formulate the problem in 
terms of the complete state variable vector 7 = (p, p, /i, a, to, X\, . . . , Xt) = (0, X). In 
this case the joint density is 

f{d t ,r t \x t ,e) = f{dt\x t ,0)f{f t \d t ,x t ,d), (19) 
and the joint likelihood function is 

T 

t=l 

The under-flow problem in evaluating (|2"U|) can now be readily overcome by the 
usual approach of working with the log-likelihood function instead (also, numerical 
integration is not required). 

By considering default and recovery processes separately and assuming a large num- 
ber of firms in the portfolio, some approximation can be justified to simplify the evalu- 
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ation of the likelihood function and its maximization procedures. This is the approach 
taken by Diillmann and Trapp (2004) and it is discussed in the next section. 

3.2 Approximate likelihood functions 

The parameter estimation procedure adopted by Diillmann and Trapp (2004) involves 
some approximations and follows two stages: in the first stage parameters for the 
default process p and p and systematic factor X are estimated; and in the second 
stage the parameters of the recovery rates are estimated. Closed-form solution can be 
found for the first stage. For the second stage, however, Diillmann and Trapp found 
that searching numerically for the maximum likelihood of the recovery model may 
provide spurious results. Thus they took a "feasible maximum likelihood" approach 
that involves two steps to estimate the recovery parameters. In the first step the 
volatility parameter a is estimated by the historical volatility a^, then in the second 
step parameters p and u are estimated conditional on o^. 

Remark: Alternatively, as discussed in the next section, the Bayesian MCMC allows 
us to estimate all parameters jointly, covering both default and recovery processes and 
systematic latent factors. The economic capital implications of taking the simplified 
two-stage and two-step maximum likelihood approach can be quantified by comparing 
the MLE of the EC value with that of the MCMC using the full joint likelihood function 



Default process 

Given X t in time period t, the conditional default probability A t = A(X t ) is a function 
of Xt given by (fl2l) . Since the function A(X) is monotonic in X, it can be inverted 
with respect to X. From f fl2|) we obtain the inverse 



7 - y/1 - p5 t 

VP 



(21) 



where 7 = $ 1 (p) and 5 t = Q l {^t)- Given that the density of X t is standard normal, 
the change of probability measure gives the density for A t : 



/(A( | 9D) = _L exp (_£T) 



dxt 



(22) 



where Op — {p, p) is the parameter vector for default process. The density of the 
transformed conditional default probability $ _1 (A t ), is then 
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V_/^ IJ 1 I 

p V 2 P ; 

For time period t we observe default rate ip t for the random variable ty t (see Section 
3.1). In the limit J t —> oo, approaches the conditional default probability A(X t ). 
Therefore approximately we can set S t = which is equivalent to observing St. 

Then the likelihood function for default process is 

T 

W D ) = l[f(S t \e D ), (24) 
t=l 

where S = (Si,...,St) is computed from the observed data vector for default rate 
if) = (if)i, . . . , i[)t), i-e. 5j = $ _1 (^i), z = 1, . . . , T. Maximizing ( 124"]) gives the following 
MLEs for p and p: 

(25) 



2 ' 




P^^^l, (26) 

where S = Ylt=i^t/T and cr| = 5^^=i(^t ~~ $) 2 /T- For details of derivation of ( 1251) 
and (|26|) . see Diillmann and Trapp (2004). The systematic factor X t is then estimated 
using (l2Tj) with default parameters p, p replaced with MLEs as 



= ™ , (27) 

VP 



where 7 = $ 



Recovery rates 

As described in Section 3.1, the average recovery rate R t in the time period t is from 
normal distribution N({Ir, a R ) with mean fi R = p + o^fbjX t and variance = cr 2 (l — 
oj)/d t . The density of in terms of recovery parameters 6 R = (p, a, a;), is obtained 
by rewriting ( ITTjl as 



n-\a • • '• I dt{r t - p- a^ujx t ) 2 \ 
f{r t \0 R ,x t ) = J- — — -exp — — . (2f 



2na 2 (l-uj) ^\ 2a 2 (l-w) 

Given T observations of the average recovery rates r = (r\, . . . ,¥t), the likelihood 
function for the recovery process is 

T 

£R(e R ,x) = l[f(r t \0 Rj x t ). (29) 
t=i 
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Diillmann and Trapp (2004) estimate Or by MLEs via maximization of f l2"9~j) with 
respect to 6 R , where x t is replaced with x t given in (f27|) . Due to numerical difficulties 
Diillmann and Trapp (2004) took a two-step approach to maximize the above log- 
likelihood function: the parameter a is estimated first by historical volatility Oh of the 
recovery rates, and then parameters p and u are estimated conditional on 5^. 

Remark: In the ML estimation performed by Diillmann and Trapp (2004), the default 
likelihood function £d{9d) is used for the estimation of (p, p) and X t , and then the 
recovery likelihood function £r(0r, x) is used for the estimation of (p, a, u>), conditional 
on estimated (p, p) and x, and with a estimated separately from historical volatility. 
For Bayesian inference, however, the joint likelihood function ffl8|) of both default and 
recovery processes has to be used for the inference of the model parameters; or use (120]) 
if X t are treated as extra parameters. 



4 Bayesian inference of LGD model using MCMC 

In the current one-factor model, the systematic risk factor X t ,t = 1, . . . , T for the ob- 
served period is latent random variable in nature while parameters are deterministic 
but unknown. From a Bayesian point of view, both and X are considered to be 
random variables; denote 7 = (0, X). Given a prior density 71(7) and a likelihood 
^(yll) — ^V(t) of a U available data Y, the joint density of Y and 7 is 

7r(y,7) = 7r(y|7M7). (30) 

Having observed data Y, the density of 7 conditional on Y (the posterior density) is 
determined by the Bayes theorem 

*(7|y) = 7r(y|7)7r(7) (31) 
J 7r(y 1 7)^(7)^7 

The posterior can then be used for predictive inference, allowing robust analysis such 
as quantification of parameter uncertainty and model choices. There is a large number 
of useful texts on Bayesian inference; for example, see Robert and Smith (1994), Lee 
(1997), Berger (1985), Robert (2001), Winkler (2003), Gelman et al (2003), Bolstad 
(2004) and Carlin and Louis (2008). In particular, examples of applying Bayesian 
inference in operational risk and insurance modeling are found in Shevchenko and 
Temnov (2009) and Peters et al (2009a, 2009b). 

The explicit evaluation of the integration in ( )3T|) is often complex especially in high 
dimensions. The Markov chain Monte Carlo (MCMC) method provides a highly ef- 
ficient alternative to traditional techniques by sampling from the posterior indirectly 
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and performing the integration implicitly. It is especially suited to Bayesian infer- 
ence framework, as it allows a unified estimation procedure, simultaneously estimating 
parameters and latent variables, and it facilitates the quantification of parameter un- 
certainty and model risks. 

4.1 Estimates of capital with parameter uncertainty 

Bayesian methods are particularly convenient to quantify parameter uncertainty and 
its impact on capital estimate. Under the Bayesian approach, model parameters 6 
are modelled by the random variable following the posterior density 7r(0|y), where 
Y = y are available data. In this case, the full predictive density (accounting for 
parameter uncertainty) of the next time period loss L T+ i, given data Y, is 

h T+ My) = j /z T+1 (z|0)7r(0|y)d0. (32) 

Here, it is assumed that, given ©, Lt+i and Y are independent. The quantile of the 
full predictive density (|32|) . 

Q p q = inf{z : Pr[L r+1 > z\Y] < 1 - q}, (33) 

at the level q, can be used as a risk measure for capital calculations; see Shevchenko 
(2008) for example of using (|33l) for operational risk capital charge. Here, "P" in 
the upper script is used to emphasize that this is a quantile of the full predictive 
distribution. The procedure for simulating from f )32|) and calculating Q q using 

posterior samples of parameters can be described as follows. 

Algorithm 2 (Quantile of full predictive loss distribution) 

1. Draw a sample 6 for the parameters from the posterior density ii(6\y) (an efficient 
sampling technique is MCMC, and the details of which will be described in Section 
4.2). 

2. Given posterior sample 6 for the parameters, simulate loss L following steps 1 to 
4 in Algorithm 1. 

3. Repeat the above steps 1-2 to obtain N samples of L. 

4. Estimate Q q using obtained samples of L in the standard way. 

Distribution of capital estimate. Another approach under a Bayesian framework 
to account for parameter uncertainty is to consider a quantile Q q (@) of the conditional 
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loss density /(-|0): 



Q g (@) = inf{z : Pr[L T+1 > z\@] <l-q}. 



(34) 



Then, given that is distributed as 7r(0|y), one can find the associated distribution 
of Q q (@) and form a predictive interval to contain the true quantile value with some 
probability. Under this approach, one can argue that the conservative estimate of the 
capital accounting for parameter uncertainty should be based on the upper bound of 
the constructed predictive interval. The following algorithm can be used to obtain the 
posterior distribution of quantile Q q (&). 

Algorithm 3 (Distribution of quantile) 

1. Draw a sample 6 for the parameters from the posterior density 7r(0\y). This can 
be done using MCMC described in Section 4.2. 

2. Compute Q q = Q g (0) using e.g. Algorithm 1. 

3. Repeat the above steps to obtain N samples of Q q (&). 

In practice the above simulation procedure for estimating the distribution of Q q (&) 
can be time consuming, because it involves a long loop (Algorithm 1) inside the loop 
over parameter samples. However, for some cases, the inner loop (Step 2) can be 
approximated by a closed-form formula and thus making the calculation of the distri- 
bution of Q q {&) more affordable. In particular, for economic capital of fully diversified 
portfolio pip , the inner loop (step 2) for computing quantile given parameters is re- 
placed by Q™{0) = A(X)S(X), with X = $^(1 - a), A(X) given by {T2]) and S(X) 

by <m>- 

4.2 Metropolis-Hastings algorithm 

Sampling from the posterior density 7r(7|y) oc 7r(y|7)7r(7) can be achieved using 
Metropolis-Hastings algorithm first described by Hastings (1970) as a generalization of 
the Metropolis algorithm (Metropolis et al 1953). Denote the state vector 7 at step m 
as 7^ m ), then the algorithm can be described as the follows. 

Algorithm 4 (Metropolis-Hastings algorithm) 

1. Start with an arbitrary initial value 7^ for m = 0. 

2. Generate 7* from the proposal density g(7*|7 (m ' ) ). 



3. Compute acceptance probability 0(7^,7*) 



min < 1 



?r(7 (m) |y)<?(7*|7 (m) ) 



n(Y\y)gh <m) h*) 



} 
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4. Draw u from the uniform distribution [7(0,1), and let ^( m+1 ) = 7* if -u < 
a(7 (m) ,7*), otherwise 7 ( m+1 ) = 7 ( m ). 

5. Repeat Steps 2 to 4 to obtain posterior samples for state variable vector 7 (col- 
lecting after burn- in period). 

The single-component Metropolis-Hastings is often efficient in practice, where the 
state variable components (7^ , 7^ , • • • , 7rf^) are updated one by one or block by 
block. This was the framework for MCMC originally proposed by Metropolis et al. 
(1953), and is adapted in this study. Other alternative MCMC methods also exist, 
e.g. the univariate slice sampler utilized by Peters et al (2009b) for estimating model 
parameters and latent factors in the context of operational risk model. 

Prior distributions 

In all the MCMC simulation runs, we assume a uniform prior for all model parameter 
6, and the standard normal distribution as the prior for latent variables Xi, . . . ,Xt- 
The only subjective judgement we bring to the prior is the lower and upper bounds of 
the parameter values. The range of the parameter value should be sufficiently large to 
allow the posterior to be implied mainly by the observed data. Table 1 shows the lower 
and upper bounds for the uniform distribution of each model parameter. As shown in 
the table, all parameters have lower and upper bounds either corresponding to the full 
support of the parameter domain or covering a wide range, and we checked that an 
increase in bounds did not lead to material difference in results. 

The starting value of Markov chain for the fcth component is set to a uniform ran- 
dom number drawn independently from the support (a K ,6 K ). In the single-component 
Metropolis-Hastings algorithm, we adopt a truncated Gaussian distribution as the sym- 
metric random walk proposal density. In addition, the Gaussian density was truncated 
below dfc and above b k to ensure each proposal is drawn within the support for the 
parameters. Specifically, for the k th component at chain step m, the proposal density 
is 



3*. /)("») rr RW\ 

q k (9*\9D = , T ' k ' k ' , , , (35) 



TP (h ■ fl^ rrRW\ J? ( n . fj( m ) ,-r.RW\ 

±<N{0k,t> k ,a k )-i<N{a k: V k ,a k 



where /jv(-; #£ m \ and F N {-;6^ k m \a k w ) are the Gaussian density and distribution 

functions respectively, with 6^ as the mean and a k w as the standard deviation. For 
each component the mean of the Gaussian density was set to the current state and the 
variance was pre-tuned and adjusted so as to allow the acceptance rate to stay at or close 
to the optimal level. For d-dimensional target distributions with i.i.d. components, the 
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asymptotic optimal acceptance rate was found to be 0.234 (Gelman et al 1997, Roberts 
and Rosenthal 2001). In pre-tuning the variances for all the components we set 0.234 
as the target acceptance rate. 

The MCMC run consists of three stages: in the first stage we tune and adjust the pro- 
posal standard deviation aj} w to achieve optimal acceptance rate for each component; 
the second stage is the "burn-in" stage and samples from this period are discarded; 
and the last stage is the posterior sampling stage where the Markov chian is considered 
to have converged to the stationary target distribution. Unless stated otherwise, the 
MCMC was performed for a length of Nb = 20, 000 as the "burn-in" period, we then 
let the chain run for an additional length of iV = 100, 000 to generate the posterior 
samples. Each iteration contains a complete update of all components. 

5 MCMC results 

5.1 Confirmation of maximum likelihood results 

We first validated our MCMC algorithm by simulating data from default and recovery 
models assuming some realistic parameter values, and performing MCMC on the sim- 
ulated data. The posterior mode should approach the assumed parameter values used 
in the simulation when the sample size increases. Having satisfied with this validation, 
we then proceed to confirm the MLE results of Dullmann and Trapp (2004), using the 
same real data and approximate likelihood functions. 

The data source for default frequencies and recovery rates is Standard&Poor's Credit 
Pro database. The actual data used by Dullmann and Trapp (2004) were kindly sup- 
plied to us by these authors. This dataset contains default rate and recovery rate 
observations for 18 years (T = 18), from January 1982 to December 1999. The recov- 
ery rates are measured either by market prices at default or prices at emergence from 
default. Dullmann and Trapp (2004) found the estimates of the expected recovery rates 
are 9% - 26% higher for prices at emergence than for market prices at default. In this 
study we use the definition with market prices at default. A very detailed description 
of the raw data and their pre-processing was given by Dullmann and Trapp (2004). 

As discussed previously, the maximum likelihood procedure adopted by Dullmann 
and Trapp (2004) involves separate stages for the default and recovery processes, and 
an approximation in the second stage by using historical volatility. To confirm their 
MLE results with our MCMC simulations, we follow the same steps. That is, we 
first performed MCMC using the default probability likelihood (|24|) . The posterior 
mode for p and p indeed agrees with those obtained using ML. Then in the second 
stage we perform MCMC with the likelihood function for recovery ( 1291) . conditioning 
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on p, p and the same historical volatility as used by Diillmann and Trapp. Again, 
the posterior modes for /i and u agree with the MLEs. These MLE results and the 
corresponding values for the the stressed probability of default PD, the stressed LGD 
and the corresponding economic capital EC (see Section 2.3.3 for definitions of these 
quantities), with X = (0.001), are shown in Table 2. 

Thereafter, unless otherwise stated, MCMC results shown are for the full joint like- 
lihood function fTSOj) . i.e. without the approximations discussed here and in Section 
3.2. 

5.2 Posterior distributions 

After validating our numerical algorithm, the full MCMC simulation was run, using 
likelihood function fTSOj) for default probability parameters (p, p) and recovery parame- 
ters (ji, cr, u), treating latent variables X = (X 1; . . . , X T ) as parameters and using data 
described above. 

The posterior sample paths (after the burn-in) for parameters (p, p, p, u, a) are 
shown in Figure 1. All paths reveal well- mixed MCMC samples indicative of sta- 
tionary distributions, as should be expected from a convergent MCMC simulation. 
Figure 2 shows the posterior density functions estimated from the posterior samples 
for parameters (p, p, p,u,a) and one of the 18 latent variables X 10 . Clearly, the den- 
sities show some positive skewness for all parameters and negative skewness for X w . 
Table 3 shows the summary statistics of all five parameters computed from the pos- 
terior samples - the mode, mean, standard deviation (stdev), skewness and kurtosis. 
To quantify uncertainty in a simple manner, the coefficient of variation (CV), defined 
as the ratio of standard deviation to the mean, is also shown in Table 3. Consistent 
with Figure 2, we see significant positive skewness in most parameters. In addition, 
the kurtosis values of some parameters are significantly higher than 3 (kurtosis of a 
normal density). 

5.3 Comparison with MLE results 

As discussed in Section 3.2, the parameter estimation procedure adopted by Diillmann 
and Trapp (2004) involves some approximations in the likelihood functions and follows 
a "feasible maximum likelihood" approach to overcome some difficulties in fitting model 
to the observed data. One of the key simplification of the two-step MLE procedure is 
that, under the assumption of a large homogeneous portfolio, the default parameters 
as well as the common factor X t are inferred from default data only, rather than jointly 
determined with recovery data. Another simplification is using a historical volatility 
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for recovery parameter a. It is interesting to see the difference between the MLE results 
with those simplifications and the full MCMC results with no approximations in the 
formulation, i.e. with a proper joint likelihood function. 

MCMC gives samples from posterior distributions for parameters 6 and latent vari- 
able (common factor) X. A meaningful comparison is between the summary statistics 
of posterior distributions with the point estimates of MLEs. Comparison between Ta- 
ble 2 and 3 shows the approximate MLEs for parameters p, p, p are within one standard 
deviation of the posterior mean of the full MCMC results, while there is a significant 
difference between MLEs and the posterior mean for the recovery parameters u and 
a. In particular, the posterior mean of recovery volatility a is more than five times 
higher than the historical value used in the approximate "feasible maximum likelihood" 
approach for the MLEs. 

A very large difference in model parameters does not always imply a large difference 
in model predictions. The comparison for the 18 latent variables (for 18 years of data) 
X t ,t = 1, . . . , 18 is shown in Figure |3j The X t implied by MLE parameter values 
using f l2~T]) agrees well with the posterior sample mean of MCMC - only one point is 
more than one standard deviation away from the posterior mean. The predictions 
on PD, LGD and EC under stressed conditions are shown in Table 2 for MLE and 
Table 4 for MCMC. The quantiles in Table 4 were obtained from Algorithm 3. For a 
comparison between point estimates, the point estimates for function PD, LGD and EC 
using the posterior mean (or MMSE, the minimum mean square estimator) Q MMSE = 
E[0\Y] of the parameters (p, p, /i,u,a) (instead of MLEs) are PD(0 MMSE ) = 0.067, 
LGD(9 MMSE ) = 0.798 and QS° 999 (6> MMSE ) = 0.0534. 

The posterior density of EC is shown in Figure 4. Evidently the distribution is 
positively skewed. Comparison of Table 2 and 4 shows that for the EC as defined 
in ( TTTT) . the MLE point estimate is 75% lower than the posterior mean from the full 
MCMC, 56% lower than the posterior median and more than 100% lower than the 0.75 
quantile. The underestimation of EC by the approximate MLE in comparison with 
full MCMC Bayesian inference is quite significant, and this is the consequence of MLE 
ignoring parameter uncertainty as well as making simplifications in the model. 

5.4 Capital estimate using full predictive loss distribution 

The 0.999 quantiles Qq 999 of the full predictive loss distribution for several portfolios 
are shown in Table 5. Here, equal weights Wj = 1/ J, j = 1, . . . , J are assumed and the 
portfolios differ only in the number of borrowers J. The highest number J = 5350 is 
the actual number of firms in the last year of the 18 years with observed data. The 
qauntiles Qq 999 in Table 5 were computed using Algorithm 2. For J = oo, instead 
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of using Step 2 of Algorithm 2, the loss L is computed using formula (ITUj) for the 
limiting case of a large portfolio. Clearly from Table 5, the full predictive quantile 
Q^999 decreases with the number of firms J, reaching a limiting value at J = oo. 
The smaller quantile for the loss distribution of a larger portfolio is a diversification 
effect. Higher risks are assumed if the loan portfolio has a small number of borrowers, 
assuming homogeneous default and recovery characteristics for each borrower. Also, 
results in Table 5 indicate that approximating the portfolio loss L in (1) by the limiting 
case formula fflUj) is valid only if the number of borrowers is in the order of thousands. 
For instance at J = 500 the full predictive quantile Q^ggg is still 10% larger than the 
limiting case for J = oo. 

The posterior density of the full predictive distribution for L°° is shown in Figure 5. 
The full predictive quantile Qq 999 at J = oo is more than twice as large as the EC 00 
estimated by the approximate MLE (shown in Table 2), and it is also 20% larger than 
the mean of Q£°(0) (shown in Table 4). This illustrates that parameter uncertainty is 
very significant in determining proper capital charge in the one-factor credit risk model 
studied here. 

6 Conclusion 

We have investigated the impact of parameter uncertainty on economic capital estimate 
in the well-known extended LGD model with systematic dependence between default 
and recovery. Uncertainty is quantified by using summary statistics of the full posterior 
distributions of model parameters obtained through Bayesian inference using Markov 
chain Monte Carlo. Posterior distribution of <5o°999(®) ( or EC) shows even the 25% 
quantile of Qo 999(©) is larger than the maximum likelihood estimated EC value; the 
mean value of the posterior distribution for EC is 75% higher than MLE estimate. In 
addition, the full predictive quantile Q^ggg for the loss distribution is more than twice 
as large as the approximate MLE estimate for the economic capital. 

The MCMC algorithm described in this work can be readily extended to deal with 
multi-factor models. Overall, it is important to account for parameter uncertainty for 
credit risk models, though it inevitably leads to unpleasant larger capital charge for 
the banks. MCMC and Bayesian approach are convenient methods to accomplish this 
task. Model considered in this study is homogeneous, but it is not difficult to apply 
the presented Bayesian MCMC for non-homogeneous case. 
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Table 1: Lower and upper bounds of model parameters. 
p p P a uj 

(0,1) (0,1) (0,1) (0.01,1.0) (0,1) 



Table 2: Maximum likelihood estimates of parameters for probability of default and 
recovery rate models and corresponding estimates of stressed PD, LGD and EC using 
approximate likelihood function ( J24l) for default and (l29l) for recovery. Recovery pa- 
rameter a is estimated by historical volatility. Here EC is estimated as Qo°999(^ MLE ) 
in (HU. 



p 


P 


P 


a co PD 


LGD 


EC 


0.0123 


0.0406 


0.438 


0.0845 0.0998 0.0476 


0.644 


0.0307 



Table 3: Summary statistics of parameters (p, p, /z, u, a) from posterior samples. Stdev 
is the standard deviation, and CV is the coefficient of variation. 





Mode 


Mean 


Stdev 


Skewness 


Kurtosis 


CV 


p 


0.0157 


0.0133 


0.0022 


0.951 


4.86 


0.168 


P 


0.143 


0.0623 


0.0239 


1.07 


4.74 


0.376 


p 


0.471 


0.456 


0.027 


0.221 


3.64 


0.058 


UJ 


0.060 


0.032 


0.023 


1.72 


8.32 


0.711 


a 


0.448 


0.457 


0.085 


0.912 


4.50 


0.183 
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Table 4: Summary statistics of stressed (assuming systematic factor x t = (0.001)) 
PD, LGD and EC posterior samples. SEC(%) is the relative difference in percentage 
between each quantile value for the distribution of Q^(G) calculated using Algorithm 
3 and the EC value estimated by MLE, Qq^qqq(0 mle ) (shown in Table 2). 



item 


Mean 


Stdev 


0.25Q 


0.5Q 


0.75Q 


CV 


PD 


0.0682 


0.0236 


0.0513 


0.0629 


0.0800 


0.346 


LGD 


0.776 


0.0769 


0.721 


0.768 


0.821 


0.099 


^0°999(®) 


0.054 


0.023 


0.0378 


0.0482 


0.0645 


0.423 


SEC(%) 


75.8% 


N/A 


23.3% 


56.9% 


110.0% 


N/A 



Table 5: Full predictive quantile Q 0.999 f° r various portfolios (different number of bor- 
rowers) using Algorithm 2. 

J = 50 J = 500 J = 5350 J = oo 

0.1044 0.0732 0.0674 0.0666 
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Figure 1: Markov chain sample paths for parameters (p, p, /i,u>,a) 
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Figure 3: Posterior mean for X t ,t = 1, . . . , 18 (solid line) in comparison with point 
estimates of X t predicted by the approximate MLEs (dots), corresponding to the 18 
years of data. Error bars correspond to posterior standard deviation of X t . 
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Figure 4: Posterior density of EC°° = Qo°999(®) computed from MCMC samples using 
Algorithm 3. 
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Figure 5: The full predictive density of L^ +1 , computed from MCMC samples using 
Algorithm 2. Dashed line indicates Q^ 999 - 
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